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We present results of our pulsar population synthesis of normal pulsars from the Galactic disk 
using our previously developed Monte-Carlo code. From our studies of observed radio pulsars that 
have clearly identifiable core and cone components, in which we fit the polarization sweep as well 
as the pulse profiles in order to constrain the viewing geometry, we develop a model describing 
the luminosity and ratio of radio core-to-cone peak fluxes. In this model, short period pulsars are 
more cone-dominated. We explore models of neutron star evolution with and without magnetic 
field decay, and with different initial period distributions. We present preliminary results including 
simulated population statistics that are compared with the observed radio pulsar population. The 
evolved neutron star populations resulting from this simulation can be used to model distributions 
of 7-ray pulsars for comparison to Fermi results. 



I. INTRODUCTION 

NASA's new flagship in 7-ray astronomy, the Fermi 
Gamma- Ray Space Telescope (Fermi) launched in 
June 2008, has already dramatically improved our un- 
derstanding of 7-ray emission from pulsars. In its first 
year of operations, Fermi has discovered over 40 new 
pulsars above 20 MeV in photon energy, superseding 
the 7-ray pulsar database of six provided by its prede- 
cessor, the EGRET instrument aboard the Compton 
Gamma- Ray Observatory. Over a dozen of these have 
been discovered through blind searches in 7 rays, most 
of which display, as yet, no evidence of a radio signal. 
In addition, Fermi has identified nearly a dozen mil- 
lisecond pulsars that have been observed in the radio. 

With the wealth of new data from Fermi at hand, 
we have the first real opportunity in decades to finally 
understand the high-energy emission and acceleration 
in pulsar magnetospheres. Making full use of this op- 
portunity will require more detailed modeling of dif- 
ferent emission mechanisms and geometries, and com- 
paring model predictions with the well-defined trends 
in the observations should be very productive. In ad- 
dition to better defining pulsar 7-ray emission, it will 
be possible also to constrain the radio emission lo- 
cation and geometry. The large collection of radio- 
quiet / radio- weak 7-ray pulsars are defining the view- 
ing angles where we are just crossing the outer edge 
of the radio beam or missing it altogether. Modeling 
both the radio beam and 7-ray beam geometry to- 
gether (as presented in companion paper Pierbattista 
et al. and was recently done for millisecond pulsars 
by Venter, Harding & Guillemot Q.) will be able to 



provide meaningful constraints of beam pattern and 
emission altitude range. 



II. MAGNETIC FIELD AND PERIOD BIRTH 
DISTRIBUTIONS 

A. Case A Magnetic field decay 

We simulate the population of normal radio pul- 
sars using two different sets of assumptions. Case A 
assumes that the pulsar spin-down can be mimicked 
through the incorporation of the decay of the mag- 
netic field with a decay constant of 2.8 Myr following 
the study of Gonthier et al. Q. While we do not ad- 
vocate that this result provides clear evidence for field 
decay, we find that this method allows one to incor- 
porate an alternative to the standard vacuum, dipole 
spin-down, as for example in the work of Contopoulos 
& Spitkovsky Q . The initial magnetic field is assumed 
to be described by two log-normal distributions given 
by the expression 



P (log B ) = ^e"( log B °- log s * 



(1) 



with the following parameters indicated in Table fl] 
The initial period distribution for this case A is given 
by a Gaussian with a mean P a and width ap o in 



P(P )oc e -( p °- p °) 2 /° 
P = 300ms 
o p = 300ms 



(2) 
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B. Case B No field decay 

Since the short decay constant of 2.8 Myr assumed 
for case A is physically unrealistic for magnetic field 
decay, we develop a no-field decay model exploring 
a radio luminosity law that is proportional to the 
square root of the spin-down power as suggested by 
the study of Faucher-Gigure & Kaspi 0]. However in 
their study, they did not include the geometry of the 
radio beam but rather used a standard beaming model 
to account for the average beam characteristics. With 
our set of assumptions defining the radio beam geome- 
try and luminosity, we are unable to reproduce the ob- 
served P — P distribution. In order to find reasonable 

I 



with the following parameters 



Miogs = 12.9 

fi Po = 200 ms (4) 
a P — 100 ms 
p = -0.6 



In Figure [TJ we display the effect of this correlated 
distribution on the present-day distribution of pulsars 
in the P — P diagram. The distribution on the left 
shows the initial distribution (red) and the present- 
day distribution (blue) of simulated pulsars assuming 
there is no correlation between the initial period and 
the magnetic field, while the distribution on the right 
shows the corresponding distributions with the cor- 
relation parameter of p = —0.6 as indicated above 
(Eq. [4]) . The distribution of detected pulsars resem- 
bles an upside-down pear shape. We find that it is 
not trivial to reproduce the narrowing of the shape 
of the distribution at the smallest with a paucity of 
pulsars to the left side (shorter period) of this narrow- 
ing around a period of about 0.7 seconds. Correlating 
the initial period with the magnetic field does seem 
to provide the necessary features to account for the 
observed distribution. We do realize that this corre- 
lation scheme is difficult to probe observationally. In 
both cases A and B, the minimum initial period in the 
simulation is set to 1.3 ms. 



agreement, we correlate the initial period distribution 
with the initial magnetic field, which remains constant 
during the pulsar spin-down. The distribution is de- 
scribed with a single log-normal B distribution and a 
correlated Gaussian P a distribution by the expression 



TABLE I: Magnetic field distribution parameters 



i 


A, 


log(Bi) 


Ci 


f 


0.6 


12.5 


0.65 


2 


0.3 


13.0 


0.80 



-2p{P -HP )(P -^ Po )/(a loEB a Po )]/[2(l-p 2 )]} (3) 
I 

III. RADIO BEAM GEOMETRY AND 
LUMINOSITY 

A. Core Cone Model 

As in our previous study [2J, we assume a beam 
model that is characterized by a central core beam 
and a single cone beam within the open-field volume 
with the cone beam decreasing in intensity near the 
last open field surface. The core beam is assumed to 
be a traditional Gaussian centered along the magnetic 
axis with a 1/e characteristic width [5[ . 

p core = 1°.5P- 1/2 (5) 

Following the work of Kijak & Gil @, we describe 
the hollow conal beam with a Gaussian having a char- 
acteristic beam radius given by 

o - 1° 24r 1/2 P~ 1 / 2 

Pcone — r K qJT 

TKG = 40^ 26 ^1°5 7 ^ Q - 30 

where tkg is the emission altitude in stellar radii. 
The characteristic beam radius corresponds to a conal 
beam opening angle at the 0.1% of the peak intensity 
of the observed profile. 

B. Radio Luminosity 

Following the standard candle assumption of radio 
pulsars, we use a luminosity model with the same 
form as used in the work of Arzoumanian, Chernoff 
& Cordes @ given by the expression 

iradio - ^^P Q Pf 5 mJy • kpc 2 • MHz (7) 



P(B,P )«CX P {-[(l0g B-Mlog B) 2 /ffl 2 og B + (ft-CP ) 2 /4 
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FIG. 1: P — P distribution of pulsars at birth (red) and at present (blue) illustrating the effect of correlating the magnetic 
field with the initial period using the same distributions in both cases. On the left the magnetic field distribution is the 
above distribution without the cross term while on the right they are correlated with the coefficient p = —0.6. 



TABLE II: Radio luminosity parameters 



Case 


R f 


a 


P 


A 


1.7 


-1.0 


0.35 


B 


2.0 


-1.5 


0.50 



pulsar spin-down that deviates from that of the tradi- 
tional vacuum dipole. New theoretical treatments are 
encouraging (see the contribution of Spitkovsky these 
proceedings) and may provide better insight on the 
pulsar spin-down. 



where Rf is a reduction factor that adjusts the birth 
rate for a given pair of a period and period deriva- 
tive exponents. P15 is the period derivative in units of 
10 -15 s • s _1 . The exponents aid in reproducing the 
observed P — P distribution, and the reduction fac- 
tor is adjusted to provide a neutron star birth rate 
of about 2.1 per century H. Our two cases A and 
B described above require different sets of coefficients 
shown below in Table |TTJ 

Following the suggestion [i| , we constrain the radio 
luminosity for the no-field decay case B to be pro- 
portional to the square root of the spin-down power 
to explore if the need for field decay stems from our 
choice of the functional form of the radio luminosity 
in our field-decay case A. Despite implementing the 
radio luminosity for case B, following the suggestion 
of Faucher-Gigure & Kaspi H, we find that to re- 
produce the P — P distribution we need to correlate 
the initial period with the magnetic field distribution 
as discussed in Section 2 for case B. The correlation 
scheme appears to have a somewhat similar effect as 
field decay, both of which are most likely mimicking a 



C. Core-to-cone peak flux model 

A core-to-cone peak flux ratio as used in the work 
of Arzoumanian, Chernoff & Cordes || is required to 
partition the luminosity between the core and cone 
components of the radio beam. We studied (Gonthier 
et al. [9() a group of 20 pulsars in the EPN database 
primarily from the observations of Gould & Lyne [icl | 
whose pulse profiles manifested three peaks. Using the 
Rotating Vector Model, we fit the polarization sweep 
to obtain the maximum rate of change of the posi- 
tion angle at the inflection point (dtp / d<p) max . This 
parameter is equal to the ratio of the sine of the mag- 
netic inclination angle and the sine of the impact an- 
gle. Using this parameter for a particular pulsar, we 
fit the profile shape allowing for the determination of 
the viewing geometry and thus obtain the ratio of the 
core-to-cone peak fluxes for those 20 pulsars within 
the assumed beam geometry and luminosity model. 
We find that we can describe this ratio with the fol- 
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FIG. 2: Aitoff plots (above) and P — P diagrams (below) for detected (left), simulated case A (center) and case B 
(right). The brown dashed lines represent lines of constant magnetic field (indicated) assuming dipole spin-down. The 
pink and green dotted curves are the curvature radiation and non-resonance inverse Compton scattering pair death lines, 
respectively. 



lowing broken power law. 



^pcak 



25 pi.3 I/ -o.9 ) for p < Js 

4 p_1.8^-0.9 for p > o J S 



(8) 



where v is the observing frequency in MHz. In this 
model, short period pulsars are much less core dom- 
inated than in the model used in the study of Ar- 
zoumanian, Chernoff & Cordes Q. In fact, when the 
period is less than 50 ms, the profile is cone domi- 
nated. This model described in Eqn. [8] appears some- 
what consistent with several recent findings. Crawford 
et al. [DIG3 studied a number of young pulsars, find- 
ing that the profiles have large linear polarization and 
little circular polarization, suggesting that the emis- 
sion is cone dominated. More recently Johnston & 
Weisberg fl3| have a similar conclusion from polariza- 
tion studies of 14 young pulsars, and Weltevrede & 
Johnston [3| found that the fractional polarization 
of pulsars is high E for pulsars with large spin-down 
powers and low E for low with the transition taking 
place near the curvature radiation pair death line. 



IV. RESULTS 

In Figure we present results of our population 
simulation for the two cases A and B, and the compar- 
ison with 1380 detected pulsars by a group of ten radio 



surveys including the Parkes Multibeam and Swin- 
burne Intermediate Latitude surveys. Above in Fig- 
ure[2]are the Aitoff plots of the detected pulsars (left), 
case A (field decay) (middle) and case B (no field de- 
cay) (right) and corresponding P — P diagrams below. 
Comparisons of the Aitoff plots suggest that the simu- 
lation is producing too many pulsars at high latitudes 
outside of the Galactic plane especially for case B. The 
distribution of detected pulsars in the P — P diagram 
is broad at high fields or large P's and narrows con- 
siderably at low fields or small P's . This trend in the 
detected distribution is more or less reproduced in the 
simulated distributions for the field decay case A and 
the no-field decay case B. In order to reproduce the 
narrowing of the distribution at low fields for case B, 
we found the need to correlate the initial period with 
the magnetic field as discussed in Section 2. With- 
out this correlation, the distribution at low fields is 
too broad with many more pulsars with small P's and 
short periods. The comparisons are not optimal, as 
the parameter space was not fully searched, which is 
one of the drawbacks of a Monte Carlo method. There 
are several interrelated parameters that affect the dis- 
tribution and further simulation runs might lead to 
improved agreement. 

Further comparisons of single pulsar characteristics 
are made in Figures 3 and 4 below by displaying his- 
tograms of various parameters (indicated) associated 
with the detected and simulated pulsars. Directly ob- 
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FIG. 3: Histograms of the indicated pulsar characteristics. The solid blue histograms represent the detected radio 
pulsars while the open red histograms represent the simulated pulsars for case A assuming magnetic field decay with a 
time constant of 2.8 Myr. 



served parameters are the period, period derivative, 
dispersion measure (DM), fluxes at 400 MHz (S400) 
and 1400 MHz (S1400). The distance of the detected 
pulsar is the best-estimated distance from the ATNF 
catalog, which if obtained from the DM, the NE2001 
electron density model of Cordes & Lazio [lj| is used 
to obtain the pulsar distance. We use this model to 
obtain the pulsar DM from the distance and Galactic 
location of the simulated pulsar. Having the distance, 
we obtain the pseudo-luminosity Sd 2 of the pulsar by 
multiplying the flux by the distance squared. The 
Galactic height and radial distance are obtained us- 
ing the distance and their Galactic location assuming 
that the Sun is located at 8.5 kpc from the Galactic 
center. For case A, assuming field decay, the age of 
the detected and simulated pulsars are recalculated 
assuming the decay constant of 2.8 Myr used in the 
simulation and zero initial period. Though not quite 
the actual age, both simulated and detect pulsars are 
treated in the same manner. As a result, the age dis- 



tribution for this case A reflects a shorter mean age 
than the one calculated from the characteristic age 
P/2P. For the pulsars in case B that assumes a con- 
stant magnetic field, the ages of the pulsars are the 
characteristic ages for both detected and simulated 
pulsars. While the simulation tracks the actual age of 
the pulsar for both case A and B, we cannot compare 
it to the ages of detected pulsars without knowing 
their initial periods. 

The Galactic scale height distributions of the simu- 
lated pulsars for both cases are broader than the one 
for the detected pulsars, which was suggested previ- 
ously in the Aitoff projections in Figure [2] The scale 
height in the initial birth distribution used in the sim- 
ulation of 75 pc is the same one that was used in 
previous study Q of normal pulsars and was not var- 
ied in the present simulations. The distributions of 
Galactic radial distance and the distance to Earth are 
well reproduced in case A but not quite as well in case 
B with more pulsars closer to us than those detected, 
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FIG. 4: Histograms of the indicated pulsar characteristics. The solid blue histograms represent the detected radio 
pulsars while the open red histograms represent the simulated pulsars for case B assuming constant magnetic field but 
correlating the initial magnetic field distribution with the initial period distribution as given by Eqn. [3] 



skewing the Galactic radial distance distribution to 
larger distances near 8 kpc. The bump appearing in 
the simulated dispersion measure distribution at low 
DM for case B is a reflection of the distance distribu- 
tion having too many nearby pulsars. The distance 
distribution is especially sensitive to the radio lumi- 
nosity. The larger than observed number of nearby 
pulsars suggest that the luminosity in the simulation 
is too small. However, the pulsar birth rate provides 
an additional constraint. The simulated birth rates 
were 2.1 and 0.9 per century for cases A and B, re- 
spectively. Increasing the luminosity further for case 
B would have improved the distance distribution, but 
would have lowered the birth rate to an even smaller 
value. 

The period and period derivative distributions are 
well reproduced for both cases A and B with a few 
too many short period pulsars in case A and a few 
too many long period pulsars in case B. Therefore the 
derived magnetic field and age distributions are well 



reproduced. The overall distributions of the directly 
measured and derived characteristics of the pulsars 
are better reproduce by case A than by case B. 



V. CONCLUSION 

From our simulations, using two sets of general as- 
sumptions, we find that the standard vacuum dipole 
spin-down of a pulsar is not able to reproduce the 
distribution of detected pulsars in the diagram, and 
either magnetic field decay (case A) or correlating the 
initial period with the magnetic field with no decay 
(case B) is necessary to qualitatively reproduce the 
upside-down pear-shape distribution of the detected 
pulsars. This result suggests that an alternative spin- 
down model maybe required. The evolved groups of 
neutron stars for both of these cases provide popula- 
tions of 7-ray pulsars that are compared to the char- 
acteristics of Fermi pulsars (see the contribution by 
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Pierbattista, Grenier, Harding & Gonthier, these pro- 
ceedings) . 
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